#!/usr/bin/env python
# coding: utf-8

# In[1]:


import numpy as np
import pandas as pd
import matplotlib.pyplot as plt
import seaborn as sns
import math
from scipy.stats import skew


# In[2]:


df = pd.read_csv("D:/research/recordbreak/record45/indus_record45_dec.csv")


# In[3]:


# Assuming 'year' is the name of the column containing the years
df['year'] = pd.to_datetime(df['year'], format='%Y').dt.year
df.set_index('year', inplace=True)


# In[4]:


print(df)


# In[45]:


dr_theo = df[df.index <= 2099][['dr_theo']]
print(dr_theo)


# In[46]:


fr_theo = df[df.index <= 2099][['fr_theo']]
print(fr_theo)


# In[7]:


#dr_obs_fact = df[df.index <= 2021][['dr_obs_fact']]
#print(dr_obs_fact)


# In[8]:


#fr_obs_fact = df[df.index <= 2021][['fr_obs_fact']]
#print(fr_obs_fact)


# # Decreasing

# In[9]:


dr_obs_dec = df[df.index <= 2022][['dr_obs']]
print(dr_obs_dec)


# In[10]:


fr_obs_dec = df[df.index <= 2022][['fr_obs']]
print(fr_obs_dec)


# In[12]:


ratio_obs_dec = df[df.index <= 2022][['ratio_obs_dec']]
print(ratio_obs_dec)


# In[13]:


dr_boot_dec = [f"dr_dec{i}" for i in range(1, 101)]

# Select the specified columns and store them in a new DataFrame
dr_boot_dec = df.loc[df.index <= 2099, dr_boot_dec]
print(dr_boot_dec)


# In[14]:


fr_boot_dec = [f"fr_dec{i}" for i in range(1, 101)]

# Select the specified columns and store them in a new DataFrame
fr_boot_dec = df.loc[df.index <= 2099, fr_boot_dec]
print(fr_boot_dec)


# In[72]:


percentile_5_dr = np.percentile(dr_boot_dec, 5, axis=1)
percentile_95_dr = np.percentile(dr_boot_dec, 95, axis=1)
percentile_50_dr = np.percentile(dr_boot_dec, 50, axis=1)


# In[73]:


percentile_5_fr = np.percentile(fr_boot_dec, 5, axis=1)
percentile_50_fr= np.percentile(fr_boot_dec, 50, axis=1)
percentile_95_fr = np.percentile(fr_boot_dec, 95, axis=1)


# In[74]:


# Calculate the ratio between corresponding percentiles of two datasets
ratio_5 = np.divide(percentile_5_fr, percentile_95_dr)
ratio_50 = np.divide(percentile_50_fr, percentile_50_dr)
ratio_95 = np.divide(percentile_95_fr, percentile_5_dr)


# In[75]:


ratio_obs_dec = df[df.index <= 2022][['ratio_obs_dec']]
print(ratio_obs_dec)


# # RCP 8.5

# In[76]:


df1 =pd.read_csv("D:/research/recordbreak/record85/indus_record85_dec.csv")


# In[77]:


# Assuming 'year' is the name of the column containing the years
df1['year'] = pd.to_datetime(df1['year'], format='%Y').dt.year
df1.set_index('year', inplace=True)
print(df1)


# In[78]:


dr85_boot_dec = [f"dr_dec{i}" for i in range(1, 101)]

# Select the specified columns and store them in a new DataFrame
dr85_boot_dec = df1.loc[df1.index <= 2099, dr85_boot_dec]
print(dr85_boot_dec)


# In[79]:


fr85_boot_dec = [f"fr_dec{i}" for i in range(1, 101)]

# Select the specified columns and store them in a new DataFrame
fr85_boot_dec = df1.loc[df1.index <= 2099, fr85_boot_dec]
print(fr85_boot_dec)


# In[80]:


percentile_5_dr85 = np.percentile(dr85_boot_dec, 5, axis=1)
percentile_95_dr85 = np.percentile(dr85_boot_dec, 95, axis=1)
percentile_50_dr85 = np.percentile(dr85_boot_dec, 50, axis=1)


# In[81]:


percentile_5_fr85 = np.percentile(fr85_boot_dec, 5, axis=1)
percentile_95_fr85 = np.percentile(fr85_boot_dec, 95, axis=1)
percentile_50_fr85 = np.percentile(fr85_boot_dec, 50, axis=1)


# In[ ]:


ratio85_5 = np.divide(percentile_5_fr85, percentile_95_dr85)
ratio85_50 = np.divide(percentile_50_fr85, percentile_50_dr85)
ratio85_95 = np.divide(percentile_95_fr85, percentile_5_dr85)


# # Combined

# In[47]:


plt.figure(figsize=(16, 9))
# Plotting the drought record-breaking events decay rate
plt.plot(dr_obs_dec.index, dr_obs_dec.values, label='Observed', color='black', alpha=0.8, linewidth=4)

# plt.plot(dr85_m_dec.index, dr85_m_dec, label='DNN', color='purple', linewidth=3) #, linestyle='dotted'
#for col in dr_boot_dec.columns:
    #plt.plot(dr_boot_dec.index, dr_boot_dec[col], color='teal', alpha=0.03)
# Plot 5th and 95th percentiles
plt.plot(dr_boot_dec.index, percentile_50_dr, color='teal', linewidth=3, alpha=0.9, label='RCP 4.5')
#plt.plot(dr_boot_dec.index, percentile_5_dr, color='teal', linestyle='dashed', alpha=0.6, linewidth=1.5)
#plt.plot(dr_boot_dec.index, percentile_95_dr, color='teal', linestyle='dashed', alpha=0.6, linewidth=1.5, label='5-95 CI')
plt.fill_between(dr_boot_dec.index, percentile_5_dr, percentile_95_dr, color='teal', alpha=0.2)

#for col in dr85_boot_dec.columns:
    #plt.plot(dr85_boot_dec.index, dr85_boot_dec[col], color='red', alpha=0.03)
# Plot 5th and 95th percentiles
plt.plot(dr85_boot_dec.index, percentile_50_dr85, color='red', linewidth=3, alpha=0.8, label='RCP8.5')
#plt.plot(dr85_boot_dec.index, percentile_5_dr85, color='red', linestyle='dashed', alpha=0.6, linewidth=1.5)
#plt.plot(dr85_boot_dec.index, percentile_95_dr85, color='red', linestyle='dashed', alpha=0.6, linewidth=1.5, label='5-95 CI')
plt.fill_between(dr_boot_dec.index, percentile_5_dr85, percentile_95_dr85, color='red', alpha=0.15)

#Theoratical
plt.plot(dr_theo, label='Theoretical (1/n)', color='blue', linewidth=3)
# Adding labels and title
plt.xlabel('Water Year(Apr-Mar)', fontsize=22,fontweight='bold') #fontweight='bold'
plt.ylabel('Probabaility (Log)', fontsize=22,fontweight='bold')
#plt.title('Indus drying Annual-Min Record-breaking', fontsize=20, fontweight='bold')

# Increase axis ticks size
plt.xticks(fontsize=20,fontweight='bold')
plt.yticks(fontsize=18,fontweight='bold')
# Increase the size of x and y tick lines
plt.gca().tick_params(axis='x', which='major', length=8, width=2)  # Set the desired length and width
plt.gca().tick_params(axis='y', which='major', length=8, width=2)

plt.yscale('log')
# plt.xticks(dr85_obs.index[::4])
# Set custom y-tick labels for logarithmic scale
yticks_log = [1, 0.5,0.2,0.1,0.05,0.02,0.01]  # Adjust these values based on your desired logarithmic scale
yticklabels = [str(val) for val in yticks_log]  # Convert values to strings for labeling
plt.yticks(yticks_log, yticklabels)
# Adding a legend
plt.legend( prop={'size': 14, 'weight': 'bold'})

numbers_descending = list(range(138, 0, -1))  # Generating a reversed sequence from 138 to 1
ax = plt.gca().twinx()
ax.plot(dr_theo.index, numbers_descending, label='Theoretical (n)', color='green', alpha=0.001, linewidth=2)
ax.invert_yaxis() 
ax.set_ylabel('Return Period (Log)', fontsize=22,fontweight='bold')
ax.yaxis.get_label().set_rotation(270)
ax.yaxis.labelpad = 20

ax.tick_params(axis='y', labelsize=18)
for label in ax.yaxis.get_ticklabels():
    label.set_fontweight('bold')
    

#boundry thickness
ax = plt.gca()
for spine in ax.spines.values():
    spine.set_linewidth(2)
plt.yscale('log')

yticks_log_secondary = [1, 2,5,10, 20, 50, 100, 150]  # Adjust these values based on your desired logarithmic scale
yticklabels_secondary = [str(val) for val in yticks_log_secondary]  # Convert values to strings for labeling
ax.set_yticks(yticks_log_secondary)
ax.set_yticklabels(yticklabels_secondary)
ax = plt.gca()
ax.tick_params(axis='x', which='major', length=8, width=2)  # Set the desired length, e.g., 10
ax.tick_params(axis='y', which='major', length=8, width=2)  # Set the desired length, e.g., 10

# Display the plot
plt.show()


# In[68]:


plt.figure(figsize=(20, 10))
#plt.title('Indus drying Annual-Min Record-breaking', fontsize=20, fontweight='bold')
# Plotting the drought record-breaking events decay rate
plt.plot(dr_obs_dec.index, dr_obs_dec.values, label='Observed', color='black', alpha=0.8, linewidth=5)

plt.plot(dr_theo, label='Stationary (1/n)', color='blue', linewidth=4)

# Plot 5th and 95th percentiles
plt.plot(dr_boot_dec.index, percentile_50_dr, color='teal', linewidth=4, alpha=0.9, label='RCP 4.5')
plt.fill_between(dr_boot_dec.index, percentile_5_dr, percentile_95_dr, color='teal', alpha=0.2)



plt.plot(dr85_boot_dec.index, percentile_50_dr85, color='red', linewidth=4, alpha=0.8, label='RCP 8.5')
plt.fill_between(dr_boot_dec.index, percentile_5_dr85, percentile_95_dr85, color='red', alpha=0.15)

#Theoratical

# Adding labels and title
#plt.xlabel('Water Year(Apr-Mar)', fontsize=22) #fontweight='bold'
#plt.ylabel('Probabaility (Log)', fontsize=22)


# Increase axis ticks size
plt.xticks(fontsize=30) #,fontweight='bold'
plt.yticks(fontsize=30)
# Increase the size of x and y tick lines
plt.gca().tick_params(axis='x', which='major', length=16, width=4)  # Set the desired length and width
plt.gca().tick_params(axis='y', which='major', length=16, width=4)

plt.yscale('log')
# plt.xticks(dr85_obs.index[::4])
# Set custom y-tick labels for logarithmic scale
yticks_log = [1, 0.5,0.2,0.1,0.05,0.02,0.01]  # Adjust these values based on your desired logarithmic scale
yticklabels = [str(val) for val in yticks_log]  # Convert values to strings for labeling
plt.yticks(yticks_log, yticklabels)
# Adding a legend
legend = plt.legend(loc='upper right', prop={'size': 40}, bbox_to_anchor=(0.89, 0.88), borderaxespad=0., ncol=2)
legend.get_frame().set_linewidth(4)
legend.get_frame().set_edgecolor('black')

numbers_descending = list(range(138, 0, -1))  # Generating a reversed sequence from 138 to 1
ax = plt.gca().twinx()
ax.plot(dr_theo.index, numbers_descending, label='Theoretical (n)', color='green', alpha=0.001, linewidth=2)
ax.invert_yaxis() 
#ax.set_ylabel('Return Period (Log)', fontsize=22,fontweight='bold')
ax.yaxis.get_label().set_rotation(270)
ax.yaxis.labelpad = 20

ax.tick_params(axis='y', labelsize=30)
#for label in ax.yaxis.get_ticklabels():
    #label.set_fontweight('bold')
    

#boundry thickness
ax = plt.gca()
for spine in ax.spines.values():
    spine.set_linewidth(4)
plt.yscale('log')

yticks_log_secondary = [1, 2,5,10, 20, 50, 100, 150]  # Adjust these values based on your desired logarithmic scale
yticklabels_secondary = [str(val) for val in yticks_log_secondary]  # Convert values to strings for labeling
ax.set_yticks(yticks_log_secondary)
ax.set_yticklabels(yticklabels_secondary)
ax = plt.gca()
ax.tick_params(axis='x', which='major', length=16, width=4)  # Set the desired length, e.g., 10
ax.tick_params(axis='y', which='major', length=16, width=4)  # Set the desired length, e.g., 10
plt.xlim(1960, 2100)
# Display the plot
plt.show()


# In[85]:


plt.figure(figsize=(20, 10))
# Plotting the drought record-breaking events decay rate
plt.plot(dr_obs_dec.index, dr_obs_dec.values, label='Observed', color='black', alpha=0.9, linewidth=5)

# Plot 5th and 95th percentiles
plt.plot(dr_boot_dec.index, percentile_50_dr, color='teal', linewidth=4, alpha=0.9, label='RCP 4.5')
plt.fill_between(dr_boot_dec.index, percentile_5_dr, percentile_95_dr, color='teal', alpha=0.2)

plt.plot(dr85_boot_dec.index, percentile_50_dr85, color='red', linewidth=4, alpha=0.8, label='RCP8.5')
plt.fill_between(dr_boot_dec.index, percentile_5_dr85, percentile_95_dr85, color='red', alpha=0.15)

# Theoretical
plt.plot(dr_theo, label='Theo (1/n)', color='blue', linewidth=3,alpha=0.8)

# Adding labels and title
#plt.xlabel('Water Year (Apr-Mar)', fontsize=22)
#plt.ylabel('Probability (Log)', fontsize=22)

# Increase axis ticks size
plt.xticks(fontsize=30)
plt.yticks(fontsize=30)

# Increase the size of x and y tick lines
plt.gca().tick_params(axis='x', which='major', length=16, width=4)
plt.gca().tick_params(axis='y', which='major', length=16, width=4)

plt.yscale('log')

# Set custom y-tick labels for logarithmic scale
yticks_log = [1, 0.5, 0.2, 0.1, 0.05, 0.02, 0.01]
yticklabels = [str(val) for val in yticks_log]
plt.yticks(yticks_log, yticklabels)

# Adding a legend
#legend = plt.legend(loc='upper right', prop={'size': 30}, bbox_to_anchor=(0.99, 0.98), borderaxespad=0.)
#legend.get_frame().set_linewidth(4)
#legend.get_frame().set_edgecolor('black')

# Set boundary thickness
ax = plt.gca()
for spine in ax.spines.values():
    spine.set_linewidth(4)

plt.xlim(1960, 2100)

# Display the plot
plt.show()


# In[50]:


plt.figure(figsize=(16, 9))
# Plotting the flood record-breaking events decay rate
plt.plot(fr_obs_dec.index, fr_obs_dec.values, label='Observed', color='black', alpha=0.8, linewidth=4)

#for col in fr_boot_dec.columns:
     #plt.plot(fr_boot_dec.index, fr_boot_dec[col], color='teal', alpha=0.05)

# Plot 5th and 95th percentiles
plt.plot(fr_boot_dec.index, percentile_50_fr, color='teal',  alpha=0.9,linewidth=3, label='RCP 4.5')
#plt.plot(fr_boot_dec.index, percentile_5_fr, color='teal', linestyle='dashed', alpha=0.6, linewidth=1.5, label='5-95 CI')
#plt.plot(fr_boot_dec.index, percentile_95_fr, color='teal', linestyle='dashed', alpha=0.6, linewidth=1.5)

plt.fill_between(fr_boot_dec.index, percentile_5_fr, percentile_95_fr, color='teal', alpha=0.2)

# 100 predictions from bootstrapping
#for col in fr85_boot_dec.columns:
    #plt.plot(fr85_boot_dec.index, fr85_boot_dec[col], color='red', alpha=0.05)

# Plot 5th and 95th percentiles
plt.plot(fr85_boot_dec.index, percentile_50_fr85, color='red',  alpha=0.8,linewidth=3, label='RCP8.5')
#plt.plot(fr85_boot_dec.index, percentile_5_fr85, color='red', linestyle='dashed', alpha=0.6, linewidth=1.5,label='5-95 CI')
#plt.plot(fr85_boot_dec.index, percentile_95_fr85, color='red', linestyle='dashed', alpha=0.6, linewidth=1.5)
plt.fill_between(fr_boot_dec.index, percentile_5_fr85, percentile_95_fr85, color='red', alpha=0.15)

# Theoretical
plt.plot(dr_theo, label='Theoretical (1/n)', color='blue', linewidth=3)
# Adding labels and title
plt.xlabel('Water Year(Apr-Mar)', fontsize=22,fontweight='bold') #fontweight='bold'
plt.ylabel('Probabaility (Log)', fontsize=22,fontweight='bold')
#plt.title('Indus drying Annual-Max Record-breaking', fontsize=20, fontweight='bold')

# Increase axis ticks size
plt.xticks(fontsize=20,fontweight='bold')
plt.yticks(fontsize=18,fontweight='bold')
# Increase the size of x and y tick lines
plt.gca().tick_params(axis='x', which='major', length=8, width=2)  # Set the desired length and width
plt.gca().tick_params(axis='y', which='major', length=8, width=2)

plt.yscale('log')
# plt.xticks(dr85_obs.index[::4])
# Set custom y-tick labels for logarithmic scale
yticks_log = [1, 0.5,0.2,0.1,0.05,0.02,0.01]  # Adjust these values based on your desired logarithmic scale
yticklabels = [str(val) for val in yticks_log]  # Convert values to strings for labeling
plt.yticks(yticks_log, yticklabels)
# Adding a legend
plt.legend( prop={'size': 14, 'weight': 'bold'})

numbers_descending = list(range(138, 0, -1))  # Generating a reversed sequence from 138 to 1
ax = plt.gca().twinx()
ax.plot(dr_theo.index, numbers_descending, label='Theoretical (n)', color='green', alpha=0.001, linewidth=2)
ax.invert_yaxis() 
ax.set_ylabel('Return Period (Log)', fontsize=22,fontweight='bold')
ax.yaxis.get_label().set_rotation(270)
ax.yaxis.labelpad = 20

ax.tick_params(axis='y', labelsize=18)
for label in ax.yaxis.get_ticklabels():
    label.set_fontweight('bold')
    

#boundry thickness
ax = plt.gca()
for spine in ax.spines.values():
    spine.set_linewidth(2)
plt.yscale('log')

yticks_log_secondary = [1, 2,5,10, 20, 50, 100, 150]  # Adjust these values based on your desired logarithmic scale
yticklabels_secondary = [str(val) for val in yticks_log_secondary]  # Convert values to strings for labeling
ax.set_yticks(yticks_log_secondary)
ax.set_yticklabels(yticklabels_secondary)
ax = plt.gca()
ax.tick_params(axis='x', which='major', length=8, width=2)  # Set the desired length, e.g., 10
ax.tick_params(axis='y', which='major', length=8, width=2)  # Set the desired length, e.g., 10

# Display the plot
plt.show()


# In[69]:


#plt.title('Indus drying Annual-Max Record-breaking', fontsize=20, fontweight='bold')
plt.figure(figsize=(20, 10))
# Plotting the drought record-breaking events decay rate
plt.plot(fr_obs_dec.index, fr_obs_dec.values, label='Observed', color='black', alpha=0.8, linewidth=5)

# Plot 5th and 95th percentiles
plt.plot(fr_boot_dec.index, percentile_50_fr, color='teal', linewidth=4, alpha=0.9, label='RCP 4.5')
plt.fill_between(fr_boot_dec.index, percentile_5_fr, percentile_95_fr, color='teal', alpha=0.2)

plt.plot(fr85_boot_dec.index, percentile_50_fr85, color='red', linewidth=4, alpha=0.8, label='RCP8.5')
plt.fill_between(fr_boot_dec.index, percentile_5_fr85, percentile_95_fr85, color='red', alpha=0.15)

#Theoratical
plt.plot(fr_theo, label='Theo (1/n)', color='blue', linewidth=4)
# Adding labels and title
#plt.xlabel('Water Year(Apr-Mar)', fontsize=22) #fontweight='bold'
#plt.ylabel('Probabaility (Log)', fontsize=22)


# Increase axis ticks size
plt.xticks(fontsize=30) #,fontweight='bold'
plt.yticks(fontsize=30)
# Increase the size of x and y tick lines
plt.gca().tick_params(axis='x', which='major', length=16, width=4)  # Set the desired length and width
plt.gca().tick_params(axis='y', which='major', length=16, width=4)

plt.yscale('log')
# plt.xticks(dr85_obs.index[::4])
# Set custom y-tick labels for logarithmic scale
yticks_log = [1, 0.5,0.2,0.1,0.05,0.02,0.01]  # Adjust these values based on your desired logarithmic scale
yticklabels = [str(val) for val in yticks_log]  # Convert values to strings for labeling
plt.yticks(yticks_log, yticklabels)
# Adding a legend
# Adding a legend
legend = plt.legend(loc='upper right', prop={'size': 30}, bbox_to_anchor=(0.99, 0.98), borderaxespad=0.)
legend.get_frame().set_linewidth(4)
legend.get_frame().set_edgecolor('black')

numbers_descending = list(range(138, 0, -1))  # Generating a reversed sequence from 138 to 1
ax = plt.gca().twinx()
ax.plot(dr_theo.index, numbers_descending, label='Theo (1/n)', color='green', alpha=0.001, linewidth=2)
ax.invert_yaxis() 
#ax.set_ylabel('Return Period (Log)', fontsize=22,fontweight='bold')
ax.yaxis.get_label().set_rotation(270)
ax.yaxis.labelpad = 20

ax.tick_params(axis='y', labelsize=30)
#for label in ax.yaxis.get_ticklabels():
    #label.set_fontweight('bold')
    

#boundry thickness
ax = plt.gca()
for spine in ax.spines.values():
    spine.set_linewidth(4)
plt.yscale('log')

yticks_log_secondary = [1, 2,5,10, 20, 50, 100, 150]  # Adjust these values based on your desired logarithmic scale
yticklabels_secondary = [str(val) for val in yticks_log_secondary]  # Convert values to strings for labeling
ax.set_yticks(yticks_log_secondary)
ax.set_yticklabels(yticklabels_secondary)
ax = plt.gca()
ax.tick_params(axis='x', which='major', length=16, width=4)  # Set the desired length, e.g., 10
ax.tick_params(axis='y', which='major', length=16, width=4)  # Set the desired length, e.g., 10
plt.xlim(1960, 2100)


# Display the plot
plt.show()


# In[83]:


plt.figure(figsize=(20, 10))
# Plotting the drought record-breaking events decay rate
plt.plot(fr_obs_dec.index, fr_obs_dec.values, label='Observed', color='black', alpha=0.9, linewidth=5)

# Plot 5th and 95th percentiles
plt.plot(fr_boot_dec.index, percentile_50_fr, color='teal', linewidth=4, alpha=0.9, label='RCP 4.5')
plt.fill_between(fr_boot_dec.index, percentile_5_fr, percentile_95_fr, color='teal', alpha=0.2)

plt.plot(fr85_boot_dec.index, percentile_50_fr85, color='red', linewidth=4, alpha=0.8, label='RCP8.5')
plt.fill_between(fr_boot_dec.index, percentile_5_fr85, percentile_95_fr85, color='red', alpha=0.15)

# Theoretical
plt.plot(fr_theo, label='Theo (1/n)', color='blue', linewidth=3,alpha=0.8)

# Adding labels and title
#plt.xlabel('Water Year (Apr-Mar)', fontsize=22)
#plt.ylabel('Probability (Log)', fontsize=22)

# Increase axis ticks size
plt.xticks(fontsize=30)
plt.yticks(fontsize=30)

# Increase the size of x and y tick lines
plt.gca().tick_params(axis='x', which='major', length=16, width=4)
plt.gca().tick_params(axis='y', which='major', length=16, width=4)

plt.yscale('log')

# Set custom y-tick labels for logarithmic scale
yticks_log = [1, 0.5, 0.2, 0.1, 0.05, 0.02, 0.01]
yticklabels = [str(val) for val in yticks_log]
plt.yticks(yticks_log, yticklabels)

# Adding a legend
#legend = plt.legend(loc='upper right', prop={'size': 30}, bbox_to_anchor=(0.99, 0.98), borderaxespad=0.)
#legend.get_frame().set_linewidth(4)
#legend.get_frame().set_edgecolor('black')

# Set boundary thickness
ax = plt.gca()
for spine in ax.spines.values():
    spine.set_linewidth(4)

plt.xlim(1960, 2100)

# Display the plot
plt.show()


# In[54]:


# Calculate the inverse of the values
inv_ratio_obs_dec = 1 / ratio_obs_dec.values
inv_ratio_50 = 1 / ratio_50
inv_ratio_5 = 1 / ratio_5
inv_ratio_95 = 1 / ratio_95

inv_ratio85_50 = 1 / ratio85_50
inv_ratio85_5 = 1 / ratio85_5
inv_ratio85_95 = 1 / ratio85_95


# In[82]:


plt.figure(figsize=(20, 10))

plt.plot( ratio_obs_dec.index,ratio_obs_dec.values, label='Observed', color='black',alpha=0.9,linewidth=5)
plt.plot(fr_boot_dec.index, inv_ratio_50, color='teal', linewidth=4, alpha=0.9, label='RCP 4.5')
plt.fill_between(fr_boot_dec.index, inv_ratio_5, inv_ratio_95, color='teal', alpha=0.2)
plt.plot(fr85_boot_dec.index, inv_ratio85_50, color='red', linewidth=4, alpha=0.8, label='RCP 8.5')
plt.fill_between(fr85_boot_dec.index, inv_ratio85_5, inv_ratio85_95, color='red', alpha=0.15)

# Theoretical
plt.axhline(y=1, color='blue', linestyle='solid',  linewidth=3,alpha=0.8) #label='Theoretical=1',

#plt.xlabel('Water Year(Apr-Mar)', fontsize=30) #fontweight='bold'
#plt.ylabel('Probabaility Ratio (Log)', fontsize=30)
#plt.title('Indus Wetting Annual-Min Record-breaking', fontsize=20, fontweight='bold')

# Increase axis ticks size
plt.xticks(fontsize=36) #,fontweight='bold'
plt.yticks(fontsize=36) #,fontweight='bold'
# Increase the size of x and y tick lines
plt.gca().tick_params(axis='x', which='major', length=18, width=4)  # Set the desired length and width
plt.gca().tick_params(axis='y', which='major', length=18, width=4)

# Set logarithmic scale for y-axis
plt.yscale('log')
# Set custom y-tick labels for logarithmic scale
yticks_log = [0.2, 0.5, 1, 2, 5, 10]  # Adjust these values based on your desired logarithmic scale
yticklabels = [str(val) for val in yticks_log]  # Convert values to strings for labeling
plt.yticks(yticks_log, yticklabels)

# Increase the thickness of the outer boundary
ax = plt.gca()
for spine in ax.spines.values():
    spine.set_linewidth(4)

# Adding a legend
#legend = plt.legend(loc='upper right', prop={'size': 30},bbox_to_anchor=(0.238, 0.98), borderaxespad=0.) #,ncol=3
#legend.get_frame().set_linewidth(4)
#legend.get_frame().set_edgecolor('black')

plt.ylim(0.19, 10)
plt.xlim(1960, 2100)
# Display the plot
plt.show()


# In[ ]:





# In[ ]:





# In[ ]:





# In[ ]:




